Susceptibility of the one-dimensional, dimerized Hubbard model 
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Abstract 



We show that the zero temperature susceptibility of the one-dimensional, 

^ ■ 

' dimerized Hubbard model at quarter-filling can be accurately determined on 

o' 

the basis of exact diagonalization of small clusters. The best procedure is 

in 



to perform a finite-size scaling of the spin velocity u a , and to calculate the 



susceptibility from the Luttinger liquid relation x = 2/7rii CT . We show that 

S' 

these results are reliable by comparing them with the analytical results that 



can be obtained in the weak and strong coupling limits. We have also used 
quantum Monte Carlo simulations to calculate the temperature dependence 



^ ■ of the susceptibility for parameters that should be relevant to the Bechgaard 



salts. This shows that, used together, these numerical techniques are able to 
give precise estimates of the low temperature susceptibility of realistic one- 
dimensional models of correlated electrons. 
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The properties of correlated electrons are fairly well understood in one dimension 
The low-energy physics of systems with low lying charge and spin excitations is well described 
by the Tomonaga-Luttinger fixed-point, while the properties of systems with a gap in the 
charge or in the spin sector can be inferred from those of the half-filled Hubbard model or 
of the Luther-Emery model respectively. When it comes to a precise understanding of the 
one-dimensional properties of real materials, like the organic conductors, this is not sufficient 
however. The ultimate goal is to understand both the high and low energy properties of 
a given material in terms of a few parameters describing the basic electronic processes, 
namely hopping integrals and Coulomb repulsions. For instance, for the Bechgaard salts, a 
good candidate to describe the electronic properties is given by the quarter-filled, dimerized 
Hubbard model defined by the Hamiltonian 

H = -t X ( C iaC i+ l, a + h.C.) - h ( C *> C *+V + h ' C + U Y n ^ n H W 

ieven,cr iodd, a i 

The parameters of this model are: i) a hopping integral ti for the short bonds; ii) a hop- 
ping integral t 2 (< ^1) for the long bonds; iii) an on-site repulsion U . In the following, 
energies will be measured in units of t\, and the basic dimensionless parameters are t 2 /^i 
for the dimerization and U/ti for the Coulomb interaction. To check whether this is a fair 
description of the low energy properties, one needs reliable estimates of the quantities that 
are accessible experimentally. For instance, the activated behaviour of the resistivity at low 
temperature in (TMTTF) 2 PF 6 is believed to be due to a small gap in the charge spectrum 
Q. To get a quantitative estimate of the charge gap in the model of Eq. (1) is a difficult 
- but not hopeless - problem that we have addressed in a previous paper || . Another use- 
ful probe of the low energy physics is provided by the static susceptibility, which has been 
measured quite extensively as a function of temperature for the Bechgaard salts The 
interpretation that has been proposed so far is based on analytical results obtained in the 
weak-coupling limit of the non-dimerized model. This is probably a good approximation 
for (TMTSF)2PF6 because the dimerization is small, but not for (TMTTF^PFg, in which 
case the ratio t 2 /ti is about 0.7, and a calculation that takes the dimerization into account 

2 



should be done before a reliable interpretation of the data can be tried. Dimerization is 
actually a very common phenomenon in one-dimensional conductors, and such a calculation 
is likely to be useful in other contexts as well. 

Surprisingly enough, it seems that there has not been any serious attempt at calculating 
the susceptibility of a dimerized model. Even for the Hubbard model, to which most of 
the effort has been devoted because it is soluble with the Bethe ansatz ||, the information 
that can be found in the litterature is not complete. The susceptibility has been calculated 
at zero temperature for any density by Shiba [0 with the Bethe ansatz equations, but its 
temperature dependence could be determined along these lines only at half-filling ||. As a 
consequence, Torrance et al || had to use approximate results obtained in the limiting cases 
of weak and strong on-site repulsion to analyze the susceptibility of TTF-TCNQ because 
the band-filling is not 1 but 0.59. A few years later, the temperature dependence of the 
susceptibility for the quarter-filled case and for U/t = 4 was determined through Monte 
Carlo simulations by Hirsch and Scalapino ||10|| , but this was just intended as an illustration, 



and the precision was not good enough at low temperature for these results to be used 
in the interpretation of experimental results. In fact, the only attempt at calculating the 
temperature dependence of the susceptibility away from half-filling is due to Bourbonnais 
|TT| . His calculation is based on g-ology. The fact that the system is not half-filled is taken 



into account by neglecting the umklapp term. While this calculation seems reliable for not 
too large values of U/t and at low temperatures, it clearly fails at higher temperature, the 
susceptibility diverging instead of vanishing as 1/T. 

In the present paper, we show that it is now possible to get accurate results for both 
the zero-temperature value and the temperature dependence of the susceptibility for one- 
dimensional models that do not have an exact solution by using standard numerical tech- 
niques together with analytical results in various limits. We have decided to concentrate 
on the model of Eq. (1) because it has a direct relevance to the Bechgaard salts, for which 
extensive experiemental data are available, but the same approach can in principle be used 
for any model, with hopefully the same success. 



Numerically, there are a priori several ways to calculate the susceptibility. One possibility 
is to go back to the original definition of the susceptibility as the derivative of the energy with 
respect to the magnetization evaluated for vanishing magnetization. This is not so accurate 
however because, with finite clusters, one has only values of the energy for discrete values 
of the total magnetization, and one first has to do a linear fit of the energy as a function of 
magnetization for small values of the magnetization. This procedure turns out to be rather 
arbitrary for the model of Eq. (1) when the repulsion is large. For one- dimensional systems, 
there is another way to calculate the susceptibility. Unless there is gap in the spin sector, 
the low energy spin excitations of electrons interacting through repulsive interactions can 
be described by a Luttinger liquid with a velocity u a , and the susceptibility can be obtained 
as x — 2/ 'iiu a |12| . One route to the susceptibility is then to determine u a numerically from 
exact diagonalization fo small clusters. There are again several possibilities to extract this 
velocity. When tijt\ 7^ 1, we have shown || that there is a gap in the charge sector of the 
model of Eq. (1). The ground-state energy is then expected to scale with the number of 
sites L as |TjJ 



where e is the energy density in the thermodynamic limit, and where the central charge c 
should be equal to 1. In the present case, this cannot be a good way to determine u„. First, 
there are logartihmic corrections to the coefficient of the 1/L 2 term which makes the estimate 
based on small clusters not very accurate. Second, and more importantly, this scaling 
behaviour is probably not satisfied fot the sizes we can study with exact diagonalization. 
The reason is quite simple: When the dimerization is very small, the gap A c in the charge 
sector is also very small. Now, the scaling of Eq. (2) can be observed only when the size 
of the systems is larger than the correlation length associated with this gap, which is given 
by £ ~ u c /A c , where u c is the slope of the dispersion of the charge excitations for energies 
larger than the gap. If the systems are small, the scaling of the ground-state will probably 
be closer to the formula for systems without a gap in the charge sector, which is similar to 




(2) 
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Eq. (2), but with u a replaced by u a + u c , and Eq. (2) would provide a totally unreliable 
estimate of u a . 

So it seems that the best way to determine u a is to look at the spectrum directly, and 
to determine the slope of the spin excitations. More precisely, the finite-size estimate of u a 
is given by 



where ko is the momentum of the ground-state, and where E(k; S; L) is the lowest energy in 
the subspace of states of momentum k and total spin S for a system of size L. This quantity 
is expected to go to the value u a of the infinite system with dominant corrections of order 
1/L. We have calculated u a (L) for systems with 8, 12 and 16 sites, and for several values of 
^2 and U/t\. This scaling was in most cases already very accurately satisfied. Small but 
not negligible deviations occured for large values of U/t\. As the form of the next-to- leading 
order corrections is not clear, we have used the systems with 12 and 16 sites to perform the 
1/L extrapolation. The results for \ are given in Fig. 1. 

We have performed several checks to convince ourselves that these estimates of the zero- 
temperature susceptibility were accurate. First, when t 2 /ti = 1, the model of Eq. (1) is 
nothing but the Hubbard model. Solving numerically the equations of the Bethe ansatz, 
Shiba 0] calculated the susceptibility for several values of the density and of U. Our results 
are in perfect agreement with the results he quoted for quarter-filling. 

In the weak-coupling limit, one can use g-ology to determine the corrections to the spin 
velocity. This leads to the following expression for the susceptibility: 



and vf = Zt^j 'Jt\ +t% for the model of Eq. (1). This is again in very good agreement 
with the slopes we have obtained numerically for small values of U/t\. 

Given the uncertainty of our numerical procedure for large values oiU/ti, it is also very 
important to check the results that we have obtained in that limit. This turns out to be 



u a (L) 



E(k + 2n/L; S = 1; L) - E(k ; S = 0;L) 
2vr/L 



(3) 




(4) 
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possible along the lines we used to calculate the 1/U corrections to the charge gap ||. When 
U jt\ is very large, the ground-state is approximately given by the product of the spinless 
fermion Fermi sea with the ground-state wave function of the spin- 1/2 Heisenberg model, 
the number of spins being equal to the number of electrons in the original model ||14|| . The 
coupling between the spins is described by an effective exchange integral given by 



where 



t 1 +t 2 2{t l - 1 2 ) 2 /2\/*i*2 \ 



m =r p <* ( 6) 

Jo v 1 — q z sin </r 



is an elliptic integral of the first kind. The susceptibility is given in terms of J e g by JT 

x = Jr (7) 



This result can be seen as a generalization of the result obtained by Klein and Seitz |T5 
for the Hubbard model. To compare with our results, we have plotted \/U as a function 
of 1/U, the result for 1/U = being deduced from the previous equations. The agreement 
is satisfactory. For large U, our numerical results for \/U are slightly scattered, but their 
are consistent with the infinite U result. We consider this agreement as a strong support in 
favour of this procedure to estimate the zero temperature susceptibility of one-dimensional 
models of correlated electrons for large values of the repulsion. 

Quite independently from these results, information on the susceptibility at moderate 
to high temperatures can be obtained from Monte Carlo simulations. Experimentally, the 
temperature dependence of the susceptibility is usually known quite accurately, and its 
interpretation is not ambiguous in the sense that it gives access to the ratio U / 1, which is a 
mesure of the size of correlations, even if one does not know the value of the hopping integrals 
accurately. This should be contrasted to the low temperature value of the susceptibility: To 
give information on the size of correlation, this quantity has to be compared with estimates 
for non interacting electrons, and, for molecular conductors, these estimates are not reliable 
if they are deduced solely from band structure calculations. 
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The world-line algorithm of Hirsch et al [IB] is known to be very convenient for one- 



dimensional models, and we have used it to determine the temperature dependence of the 
susceptibility for two sets of parameters: U/t\ — 4, t 2 /£i — 1, which should be reasonable for 
(TMTSF) 2 PF 6 , and U/h = 8, t 2 /t x = -7, which should be reasonable for (TMTTF) 2 PF 6 . 
In their original study of the extended Hubbard model, Hirsch et al worked at temperatures 
such that the finite size effects were smaller than the statistical errors. In our case, we need 
to go to temperatures that are small enough to make the link with the zero temperature 
results obtained with exact diagonalization. So we had to obtain results with statistical 
errors small enough to allow meaningful extrapolations. More precisely, we had to do three 
extrapolations. First, in the algorithm we have used, the magnetization is a conserved 
quantity. So, to get the static susceptibility, one has to extrapolate the finite q results to get 
the zero-g value of the susceptibility. Practically, this is done by fitting the results at small 
q with a parabola. Then, there is the systematic error due to the Trotter decomposition 
which goes as At 2 , and which we eliminated by fitting the results obtained for different 
numbers of temperature slices. Finally, for the lowest temperatures studied, the finite-size 
effects were not negligible, and, by analogy with the procedure used for the zero-temperature 
susceptibility, we performed a linear fit in 1/L of the results obtained for different sizes to 
go to the thermodynamic limit. The error bar depicted on the figure is the largest of the 
statistical errors obtained for the various sizes and time-slice numbers for a given value of 
the temperature. The results are shown in Fig. 3. The agreement with the zero-temperature 
result is very good. This is particularly satisfactory for the case U/ti = 8, t 2 /ti = .7, because 
this was already in the region where the finite-size corrections to u a were not purely 1/L for 
the small systems we could study with exact diagonalization. 

As far as the Bechgaard salts are concerned, we will limit ourselves to a few remarks. 
A more complete account of these results, together with the interpretation of several other 
experimental data, will be given elsewhere. Let us note for the moment that the present 
results roughly confirm the interpretation given by Wzietek et al Q for (TMTSF) 2 PF 6 : 
Our result for the temperature dependence of the susceptibility for U/t\ = 4, t 2 /ti = 1 is 



consistent with the results of Bourbonnais. For (TMTTF) 2 PF 6 , our results suggest that the 
dimerization reduces the temperature dependence of the susceptibility, so that the analysis 
of Wzietek et al, which did not take the dimerization into account, leads to an underestimate 
oiU/vF in that compound. 

In conclusion, we have shown that it is possible to obtain accurate results concerning 
the susceptibility of one-dimensional models of correlated electrons that cannot be solved 
exactly by Bethe ansatz by using numerical techniques and approximate methods in various 
limits. A finite-size scaling of the spin velocity deduced from the spectrum obtained by exact 
diagonalization provides good values of the zero-temperature susceptibility, while standard 
Monte Carlo simulations give information on the temperature dependence that seems accu- 
rate enough to allow an interpretation of the experimental data obtained on the Bechgaard 
salts. It is our hope that this work will encourage experimentalists to try to interprete the 
low energy results they can obtain for one-dimensional conductors in terms of microscopic 
models. 

We acknowledge useful discussions with L. Brossard, D. Jerome, D. Poilblanc, J. -P. 
Pouget, S. Sorella, P. Wzietek and T. Ziman. We are especially indebted to M. Dzierzawa 
for providing us with a Quantum Monte Carlo program that we could adapt to calculate 
the finite temperature susceptibility. The numerical simulations were performed on the C98 
of the IDRIS (France). 
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FIGURES 

FIG. 1. Zero temperature susceptibility as a function of U/t\ for several values of the dimer- 
ization: a) £2/^1 = 0.4, 0.5, 1 from top to bottom; b) t2/t\ = 0.1, 0.2, 0.3 from top to bottom. 

FIG. 2. x^i/U as a function of t±/U for several values of the dimerization: a) £2/^1 = 0.4, 0.5, 
1 from top to bottom; b) ^Ai = 0.1, 0.2, 0.3 from top to bottom. 

FIG. 3. Temperature dependence of the susceptibility for: a) U/t\ = 4, £2/^1 = 1; b) U/t\ = 8, 

t 2 /h = 0.7. 
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